subroutine fereg(y,x,nmax,nn,maxt,kmax,k,nt,beta)
implicit double precision (a-h,o-z)
double precision :: x(nmax,maxt,kmax),y(nmax,maxt),beta(k)
integer :: nt(nmax)
double precision :: xx(k,k),xxi(k,k),xy(k),xtil(maxt,k)

xx=0.0d0
xy=0.0d0

ic=0
do i=1,nn
 if(nt(i)>1) then
  do j=1,k
   xtil(1:nt(i),j)=x(i,1:nt(i),j)-sum(x(i,1:nt(i),j))/dble(nt(i))
  enddo
  do iw=1,nt(i)
  ic=ic+1
   do j1=1,k
    xy(j1)=xy(j1)+xtil(iw,j1)*y(i,iw)
    do j2=1,k
     xx(j1,j2)=xx(j1,j2)+xtil(iw,j1)*xtil(iw,j2)
    enddo
   enddo
  enddo
 endif
enddo

call invert(xx,xxi,k)

beta=0.0d0
do j1=1,k
 do j2=1,k
  beta(j1)=beta(j1)+xxi(j1,j2)*xy(j2)
 enddo
enddo 
  
return
end
